********************************************************************************
* Purpose: Main regressions and output
********************************************************************************

do "D:\Data\workdata\708296\Project - Mother groups and mental health\Do\global.do" 

use "$work\final_data1", clear

********************************************************************************
* Figure 2: The effect of exposure to a depressed peer on Maternal mental health care usage
********************************************************************************
*regs
eststo clear
forvalue i = 1/24 {
eststo: areg anyw`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\any`i'", replace)
}
eststo clear
forvalue i = 1/24 {
eststo: areg psykw`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\psyk`i'", replace)
}
eststo clear
forvalue i = 1/24 {
eststo: areg consw`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\cons`i'", replace)
}
forvalue i = 1/24 {
eststo: areg adm`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\adm`i'", replace)
}
* graphics
use "$work\any1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\any`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\any1.png", replace

use "$work\psyk1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\psyk`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
, $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\psyk1.png", replace


use "$work\cons1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\cons`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\cons1.png", replace

use "$work\adm1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\adm`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\adm1.png", replace

********************************************************************************
* As table - Table 3 of the paper (main result)
********************************************************************************
use "$work\final_data1", clear

preserve
eststo clear
foreach var in anyw24 consw24 psykw24 adm24  {
eststo: areg `var' peer_treat , absorb(time fe mfr_first) vce(cluster group) 
estimates store p`var'
summarize `var' if peer_treat == 0, meanonly
estadd scalar means = r(mean)
}
esttab p* using "$out\maintable.tex", fragment b(3) se(3) keep(peer_treat) stats(means N, fmt(%13.2fc %13.0fc) label("Control group mean" "Obs.")) se star(* 0.10 ** 0.05 *** 0.01) alignment(c c) label nolegend nonote replace

restore

//For fathers: Appendix Table C9
eststo clear
foreach var in any2f gpf_cons2 gpf_psyk2 postbirthf2y{
eststo: areg `var' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
estimates store t`var'
summarize `var' if peer_treat == 0, meanonly
estadd scalar means = r(mean)
}
esttab t* using "$out\maintable_f.tex", fragment b(3) se(3) keep(peer_treat) stats(means N, fmt(%13.2fc %13.0fc) label("Control group mean" "Obs.")) se star(* 0.10 ** 0.05 *** 0.01) alignment(c c) label nolegend nonote replace


********************************************************************************
* Figure C6: Effects of exposure to a depressed peer on nurse care
********************************************************************************
use "$work\final_data1", clear

eststo clear
forvalues i = 1/12 {
eststo: areg unim`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\uni`i'", replace)
}
eststo clear
forvalues i = 1/12 {
eststo: areg needm`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\need`i'", replace)
}
eststo clear
forvalues i = 1/12 {
eststo: areg phonem`i' peer_treat, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\phone`i'", replace)
}
* nurse outcome graphics
use "$work\uni1",clear
gen m = 1
forvalue v = 2/12 {
append using "$work\uni`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw12 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.2(0.1)0.2, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\uni1.png", replace

use "$work\need1",clear
gen m = 1
forvalue v = 2/12 {
append using "$work\need`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw12 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.2(0.1)0.2, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\need1.png", replace

use "$work\phone1",clear
gen m = 1
forvalue v = 2/12 {
append using "$work\phone`v'"
replace m = `v' if m == .
}
keep if parm == "peer_treat"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(vsmall) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(vsmall) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw12 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.2(0.1)0.2, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\phone1.png", replace

********************************************************************************
* Nurse outcomes as table including new assignment - Appendix Table C10
********************************************************************************

use "$work\final_data1", clear
estimates clear
foreach var in unim12 needm12 phonem12 newassignment{
eststo: areg `var' peer_treat , absorb(time fe mfr_first) vce(cluster group) 
estimates store p`var'
summarize `var' if peer_treat == 0, meanonly
estadd scalar means = r(mean)
}
esttab p* using "$out\nurse_table.tex", fragment b(3) se(3) keep(peer_treat) stats(means N, fmt(%13.2fc %13.0fc) label("Control group mean" "Obs.")) se star(* 0.10 ** 0.05 *** 0.01) alignment(c c) label nolegend nonote replace

**********************************************************************************************
* Figure C3: The effect of average screening scores of peers on Mental health care usage of high scoring mothers 
**********************************************************************************************
use "$work\final_data1", clear

keep if high == 1
eststo clear
forvalue i = 1/24 {
eststo: areg anyw`i' score_e1_mean i.score_e1, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\anyhigh`i'", replace)
}

eststo clear
forvalue i = 1/24 {
eststo: areg psykw`i' score_e1_mean i.score_e1, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\psykhigh`i'", replace)
}
eststo clear
forvalue i = 1/24 {
eststo: areg consw`i' score_e1_mean i.score_e1, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\conshigh`i'", replace)
}
forvalue i = 1/24 {
eststo: areg adm`i' score_e1_mean i.score_e1, absorb(time fe mfr_first) vce(cluster group) 
parmest, format(estimate min95 max95 %8.3f p %8.1e) saving("$work\admhigh`i'", replace)
}
*************************************
* figures
*************************************
use "$work\anyhigh1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\anyhigh`v'"
replace m = `v' if m == .
}
keep if parm == "score_e1_mean"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(tiny) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(tiny) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\anyhigh1.png", replace

use "$work\psykhigh1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\psykhigh`v'"
replace m = `v' if m == .
}
keep if parm == "score_e1_mean"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(tiny) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(tiny) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\psykhigh1.png", replace

use "$work\conshigh1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\conshigh`v'"
replace m = `v' if m == .
}
keep if parm == "score_e1_mean"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(tiny) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(tiny) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\conshigh1.png", replace

use "$work\admhigh1",clear
gen m = 1
forvalue v = 2/24 {
append using "$work\admhigh`v'"
replace m = `v' if m == .
}
keep if parm == "score_e1_mean"

replace t = round(t,0.001)
replace stderr = round(stderr,0.001)
twoway (scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(estimate) mlabp(12) mlabs(tiny) mlabg(15pt) mlabc(black)) /// 
(scatter estimate m, msize(small) msymbol(Oh) mcolor(red) alwidth(none) lcolor(black) lwidth(thick) mlabel(stderr) mlabp(12) mlabs(tiny) mlabg(6pt) mlabc(black)) ///
(line estimate m, fcolor(black%10) alwidth(none) lcolor(red) lwidth(thick)) ///
(line max95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)) ///
(line min95 m, fcolor(black%30) alwidth(none) lcolor(gs10) lwidth(medium) lpattern(dash)), $gsw2 legend(off) xtitle("Age in months") ytitle("Difference") ylab(-0.03(0.01)0.03, tlcolor(gs10) labgap(6-pt) format(%4.2f)) yline(0, lcolor(gs10) lwidth(medthick))
graph export "$out\admhigh1.png", replace
